Dynamic Stability of Nanobeams Based on the Reddy’s Beam Theory

The dynamic stability of nanobeams has been investigated by the Euler-Bernoulli and Timoshenko beam theories in the literature, but the higher-order Reddy beam theory has not been applied in the dynamic stability evaluation of nanobeams. In this work, the governing equations of the motion and dynamic stability of a nanobeam embedded in elastic medium are derived based on the nonlocal theory and the Reddy’s beam theory. The parametric studies indicate that the principal instability region (PIR) moves to a lower frequency zone when length, sectional height, nonlocal parameter, Young’s modulus and mass density of the Reddy nanobeam increase. The PIR shifts to a higher frequency zone only under increasing shear modulus. Increase in length makes the width of the PIR shrink obviously, while increase in height and Young’s modulus makes the width of the PIR enlarge. The sectional width and foundation modulus have few effects on PIR.


Introduction
For a beam under an axially acted periodic excitation, when the excitation frequency is two times the natural frequency of beam, divergent oscillation can be found in the transverse direction of the beam. This phenomenon is referred to as dynamic instability. On the other hand, nanobeam structures such as carbon nanotubes (CNTs) are seen as the most promising new material to play an important role in nanotechnology [1,2], and vibration analysis of nano-scale or nano-composite structures has been carried out extensively [3,4]. The dynamic instability of nanobeams remains an open topic in recent years, and the theory of nonlocal continuum elasticity initiated by Eringen [5,6] is usually employed to take into account quantum effects at nanoscale. The major difference between the classical and nonlocal elasticity theories is that the former assumes that stress state at a given point is uniquely affected by strain state at the same point, while the latter considers stress state at a given point as a function of strain states of all points in the continuum.
Reddy [7] applied the nonlocal theory to develop various nonlocal beam models for bending, vibration, and stability of Euler-Bernoulli, Timoshenko, Reddy, and Levinson beams. Subsequently, various mechanical behaviors of nanobeams have been studied by nonlocal beam theories. Arda and Aydogdu [8] investigated the dynamic stability problem of a nanobeam under a time-varying axial loading based on the nonlocal Euler-Bernoulli beam model, which leads to a time-dependent Mathieu-Hill equation. Li et al. [9] studied the dynamics and stability of nonlocal nanobeams by the perturbation method, indicating that nonlocal nanoscale has significant effect on instability region. Ansari et al. [10] analyzed the dynamic stability of embedded single-walled carbon nanotubes (SWCNTs) in a thermal environment based on the nonlocal Bernoulli-Euler and Timoshenko beam theories.
Wang and Li [11] presented the nonlinear primary resonance of nanobeams by nonlocal continuum theory and discussed the influence of Winkler foundation modulus and the ratio of length to diameter. Ansari et al. [12] analyzed the dynamic stability of multi-walled

Governing Equation of Reddy Beam Embedded in Elastic Matrix
The mechanical behavior of a Reddy beam embedded in an elastic medium has not yet been studied in the literature, so the governing equation of motion is derived in this section.
A simply supported Reddy-type nanobeam with a rectangular cross-section is illustrated in Figure 1. The nanobeam is embedded in an elastic matrix, and the Winkler foundation modulus is k w . The nanobeam is subjected to an axial periodic excitation F(t), where t is the time variable. The length of nanobeam is represented by l, and the width and height of cross-section are denoted by b and h, respectively. For the material properties, Young's modulus, shear modulus, and mass density are represented by E, G, and ρ, respectively.
The mechanical behavior of a Reddy beam embedded in an elastic medium has not yet been studied in the literature, so the governing equation of motion is derived in this section.
A simply supported Reddy-type nanobeam with a rectangular cross-section is illustrated in Figure 1. The nanobeam is embedded in an elastic matrix, and the Winkler foundation modulus is w k . The nanobeam is subjected to an axial periodic excitation ( ) where t is the time variable. The length of nanobeam is represented by l , and the width and height of cross-section are denoted by b and h, respectively. For the material properties, Young's modulus, shear modulus, and mass density are represented by E , G , and  , respectively. A coordinate system ( , , ) x y z is introduced, where the x -coordinate is taken along the length of the nanobeam, y-coordinate along the width of the nanobeam, and z -coordinate along the thickness (the height) of the nanobeam.
The displacement field of the refined Reddy beam can be expressed as [7] where the displacements The non-zero strains in the Reddy beam can be obtained as A coordinate system (x, y, z) is introduced, where the x-coordinate is taken along the length of the nanobeam, y-coordinate along the width of the nanobeam, and z-coordinate along the thickness (the height) of the nanobeam.
The displacement field of the refined Reddy beam can be expressed as [7] where the displacements (u 1 , u 2 , u 3 ) denote the displacements along the coordinates (x, y, z); u(x, t) and w(x, t) are the axial and transverse displacements on the mid-plane (z = 0) of the beam; ϕ(x, t) is the rotation of cross-section, and The non-zero strains in the Reddy beam can be obtained as where ε xx and ε xz denote the longitudinal normal strain and the transverse shear strain, respectively. ε 0 xx is the extensional strain k R is the bending strain of Reddy beam k R = ∂ϕ ∂x (9) and α R = −c 1 ∂ϕ ∂x It can be seen from Equation (7) that the transverse shear strain vanishes on the top and bottom faces of the beam (z = ±h/2), so it no longer needs a shear correction factor as in the Timoshenko beam theory.
The kinetic energy and strain energy of Reddy beams have the forms of [7,17] (14) where N, M, P, Q, and R are the stress resultants given in Reddy [7]. The potential energies for the axial excitation and action of the elastic medium are [15,23] where F m = k w w is the force from the elastic foundation. Then, the Hamilton's principle of axially excited Reddy beams on the Winkler-type elastic foundation can be written as where the constants m i are for a rectangular cross-section nanobeam [24,25] Equation (17) leads to the Euler-Lagrange equations of motion for the Reddy beam as where

Equation of Motion Based on the Nonlocal Theory
According to the nonlocal theory that the stress at a point in an elastic continuum is the function of strains of all points in the continuum, we have the constitutive relationship in the integral form of nonlocal stress tensor [5,26] where σ ij and ε ij are the stress and strain tensors, respectively; ϑ(|x − x |, τ) is the nonlocal modulus function, denoting the nonlocal effect at point x produced by the point at x , τ is a material constant affected by internal and external characteristics length, τ = e 0 a/l, l is the external characteristic length (material length scale parameter), a is the internal characteristic length (crystal lattice parameter), e 0 is a constant identified from atomic simulation or the dispersion curve of the Born-Karman model using crystal lattice dynamics, C ijkl is the elasticity modulus tensor of classical isotropic elasticity, and u i is the displacement vector. The integral form of the nonlocal constitutive relationship is complicated and it has the simpler differential form for unidimensional nanobeams where µ = (e 0 a) 2 is the nonlocal parameter. Substituting Equations (27) and (28) and Equations (6) and (7) into the expressions of N, M, P, Q, and R, we obtain the following force-strain relations for the nonlocal Reddy beam where (A, I, J, K) = A (1, z 2 , z 4 , z 6 )dA are the second, fourth, and sixth order moments of the area around the y-axis. Further substituting Equations (29) to (32) into Equations (22) and (23) leads to the equation of motion for the nonlocal Reddy nanobeam:

Governing Equation of Dynamic Instability
For a simply supported nanobeam, its displacements w(x, t) and ϕ(x, t) can be expressed as where n denotes the number of modes. Substituting Equations (38) and (39) into Equations (33) and (34) gives where the dot superscript indicates the derivation with respect to time.
In order to ensure that Equations (38) and (39) really satisfy Equations (33) and (34), it is necessary and sufficient that the quantity in the brackets before sin(nπx/l) and cos(nπx/l) should vanish at any time t, which can be written in matrix form is displacement vector, is mass matrix, is stiffness matrix and is geometric stiffness matrix, with When the axial excitation is static, the term of inertia force in Equation (42) disappears, and Equation (42) is reduced to the problem of static instability: where F cr is Euler critical load. By disregarding axial excitation and considering dynamic displacement vector as d = de iωt , Equation (42) degenerates into the free vibration problem where ω denotes the natural frequency of the Reddy nanobeam. When F(t) is periodically time dependent, it can be defined as where α and β denote the ratios of mean and amplitude of F(t) with respect to Euler load, and θ is the excitation frequency. Substituting Equation (58) Equation (59) is a second-order differential equation with periodic coefficients, called the Mathieu-Hill equation. The solutions of the Mathieu-Hill equation, i.e., the dynamic responses of the Reddy nanobeam, can be categorized into those that are stable and those that are unstable, and these solutions cluster together to form the stable and unstable regions, respectively. The regions of stable solutions correspond to dynamic stability state of the nanobeam, so they are called the dynamic stability regions (DSRs). Meanwhile, the regions of unstable solutions correspond to dynamic instability of the nanobeam and they are named the dynamic instability regions (DIRs). Moreover, the boundaries between DSRs and DIRs are determined by the periodic solutions of the Mathieu-Hill equation with the periods of T and 2T (T = 2π/θ). Each period T or 2T corresponds to a series of DIRs, and the first DIRs under period 2T has the largest width and is the most dangerous region, which is called the principal dynamic instability region (PIR).
The boundary of the PIR can be determined by solving the following eigenvalue problems [13,15]: Since Equation (60) is identical for all n modes, the index n can be dropped in solving Equation (60). For a deterministic nanobeams, matrices K e , K g , and M and Euler load F cr are fixed, and the corresponding excitation frequency θ can be obtained under the given α and β. Thus, boundaries of the PIR can be depicted on a three-dimensional parametric space enclosed by α, β, and θ. Excitation parameter e = β/[2(1 − α)] used in the classic Bolotin's method and nondimensionalized excitation frequency θ = θL ρ/E can be further employed for conveniently demonstrating the PIR on a two-dimensional parametric plane (θ, e) [17,19].

Results and Discussions
Material and geometric parameters of the studied nonlocal Reddy beam are listed in Table 1. The corresponding Euler load and natural frequency are 1.15 × 10 −8 N, and 6.83 × 10 12 Hz, respectively, and the boundary of the PIR is given in Figure 2. According to Figure 2, the dimensionless excitation frequency is θ = 2.26 under e = 0, i.e., the excitation frequency is θ = 13.16 × 10 12 Hz, which is about two times the natural frequency. Table 1. Parameters of the studied nanobeam.

Material Parameters
Geometric Parameter   is outside the PIR, the dynamic response is stable with time.
The mathematical formulations in Section 2 have indicated that the boundary of the PIR is controlled by the geometric and material parameters shown in Table 1. In the following, the effects of each parameter on the boundary of the PIR are investigated. For the investigated parameter, other parameters in Table 1  As shown in Figure 2, when the excitation parameter is within 1.0, the Reddy nanobeam might become dynamically unstable when the dimensionless excitation frequency is in the range of [2.165, 2.347]. If the point (θ, e) is located inside the PIR, the corresponding dynamic response takes on the currency of dynamic instability. Otherwise, when the point (θ, e) is outside the PIR, the dynamic response is stable with time.
The mathematical formulations in Section 2 have indicated that the boundary of the PIR is controlled by the geometric and material parameters shown in Table 1. In the following, the effects of each parameter on the boundary of the PIR are investigated. For the investigated parameter, other parameters in Table 1 are kept invariant.
Comparison of PIRs under various lengths of the Reddy nanobeam is shown in Figure 3. It is indicated that the PIR shifts to the zone with higher frequency of excitation when the length increases and excitation frequency is nondimensionalized by θ = θL ρ/E ( Figure 3a). In order to eliminate the influence of length during nondimensionalization, excitation frequency is rescaled by θ = θ ρ/E and the PIRs are redrawn in Figure 3b. It is demonstrated that the PIR moves to the zone with lower excitation frequency when the length increases, and the width of the PIR decreases significantly.
PIR is controlled by the geometric and material parameters shown in Table 1. In the following, the effects of each parameter on the boundary of the PIR are investigated. For the investigated parameter, other parameters in Table 1 are kept invariant.
Comparison of PIRs under various lengths of the Reddy nanobeam is shown in Figure 3. It is indicated that the PIR shifts to the zone with higher frequency of excitation when the length increases and excitation frequency is nondimensionalized by Figure 3a). In order to eliminate the influence of length during nondimensionalization, excitation frequency is rescaled by E     and the PIRs are redrawn in Figure 3b. It is demonstrated that the PIR moves to the zone with lower excitation frequency when the length increases, and the width of the PIR decreases significantly. When the height of the cross-section of the nanobeam increases and θ = θL ρ/E, PIR also shifts to a lower frequency zone, and the width of the PIR increases obviously (Figure 4a). However, sectional width has a small impact on the position and width of the PIR (Figure 4b). PIR also shifts to a lower frequency zone, and the width of the PIR increases obviously (Figure 4a). However, sectional width has a small impact on the position and width of the PIR (Figure 4b).
(a) (b)       With the increase in Young's modulus, the PIR also shifts to a lower frequency zone when excitation frequency is nondimensionalized by θ = θL ρ/E (Figure 6a). For eliminating the effect of Young's modulus in the nondimensionalization, excitation frequency is rescaled by θ = θL √ ρ in Figure 6b. It is shown that width of the PIR increases when Young's modulus is raised.  In contrary to Young's modulus, the increment of shear modulus makes the PIR move to a higher frequency zone (Figure 7, θ = θL ρ/E).    Mass density has very few effects on the PIR when excitation frequency is nondimensionalized by θ = θL ρ/E (Figure 8a), but mass density makes the PIR move to a lower frequency zone when excitation frequency is zoomed by θ = θL/ √ E (Figure 8b). As illustrated in Figure 9a, foundation modulus has limited influence on the PIR when the nanobeam is short. However, with increase in length, the effect of foundation modulus becomes more remarkable, where the PIR moves to a lower frequency zone and its width becomes larger (Figure 9b). Hence, elasticity of foundation would aggravate the dynamic stability of Reddy nanobeams. As illustrated in Figure 9a, foundation modulus has limited influence on the PIR when the nanobeam is short. However, with increase in length, the effect of foundation modulus becomes more remarkable, where the PIR moves to a lower frequency zone and its width becomes larger (Figure 9b). Hence, elasticity of foundation would aggravate the dynamic stability of Reddy nanobeams. As illustrated in Figure 9a, foundation modulus has limited influence on the PIR when the nanobeam is short. However, with increase in length, the effect of foundation modulus becomes more remarkable, where the PIR moves to a lower frequency zone and its width becomes larger (Figure 9b). Hence, elasticity of foundation would aggravate the dynamic stability of Reddy nanobeams.

Conclusions
Considering that the dynamic instability of nonlocal Reddy beams has not yet been investigated, this study carefully derives the governing equations of motion and dynamic stability for a simply supported Reddy nanobeam embedded in elastic medium, based on the nonlocal theory and Reddy's beam theory.
The study demonstrates that the increase in length of a Reddy nanobeam makes the principal instability region (PIR) move to a lower frequency zone and the width of the PIR shrink obviously. Increase in sectional height of the nanobeam also causes a shift of the PIR to a lower frequency zone, but the width of the PIR enlarges significantly. Sectional width has few effects on the PIR. Under a larger nonlocal parameter, the PIR shifts to a lower frequency zone while the width of PIR is less influenced.

Conclusions
Considering that the dynamic instability of nonlocal Reddy beams has not yet been investigated, this study carefully derives the governing equations of motion and dynamic stability for a simply supported Reddy nanobeam embedded in elastic medium, based on the nonlocal theory and Reddy's beam theory.
The study demonstrates that the increase in length of a Reddy nanobeam makes the principal instability region (PIR) move to a lower frequency zone and the width of the PIR shrink obviously. Increase in sectional height of the nanobeam also causes a shift of the PIR to a lower frequency zone, but the width of the PIR enlarges significantly. Sectional width has few effects on the PIR. Under a larger nonlocal parameter, the PIR shifts to a lower frequency zone while the width of PIR is less influenced.
When the Young's modulus of the nanobeam increases, the PIR moves to a lower frequency zone and the width of the PIR enlarges. On the contrary, with increase in shear modulus, the PIR moves to a higher frequency zone and the width of the PIR shrinks. When mass density increases, the PIR moves to a lower frequency zone. Foundation modulus also has few effects on the PIR under shorter nanobeams, but the effect becomes significant when the nanobeam becomes longer.